Diamagnetic susceptibility obtained from the six-vertex model and its implications for the 
high-temperature diamagnetic state of cuprate superconductors 
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We study the diamagnetism of the 6- vertex model with the arrows as directed bond currents. To our knowl- 
edge, this is the first study of the diamagnetism of this model. A special version of this model, called F model, 
describes the thermal disordering transition of an orbital antiferromagnet, known as d-density wave (DDW), a 
proposed state for the pseudogap phase of the high-T c cuprates. We find that the F model is strongly diamag- 
netic and the susceptibility may diverge in the high temperature critical phase with power law arrow correlations. 
These results may explain the surprising recent observation of a diverging low-field diamagnetic susceptibility 
seen in some optimally doped cuprates within the DDW model of the pseudogap phase. 

PACS numbers: 73.43.Nq,74.25.Dw,75.10.Hk 



Introduction. Experiments on the normal state properties 
of the cuprate superconductors continue to pose new theoreti- 
cal challenges. Above the superconducting transition temper- 
ature T c , the cuprates in the underdoped regime evince a d- 
wave-like gap even in the absence of superconductivity. The 
nature of the system in this pseudogap phase is believed to 
hold the key [1] to the physics of the high transition tempera- 
ture itself. A recent remarkable set of experiments yl S have 
found evidence of enhanced diamagnetism in the pseudogap 
phase above T c at a doping range near and below the optimal 
doping. In particular, these experiments have revealed that, 
near optimal doping, the low-field diamagnetic susceptibility 
X diverges above T c as an inverse power of the applied field 
H, x ~ —H~ x . Here x is a T-dependent exponent. The 
divergence of x above T c implies underlying critical correla- 
tions in an entire phase above T c (01 which is not easy to ex- 
plain by any existing theories of the pseudogap phase JHIsl]- In 
this paper, we will address this question within the framework 
of the d-density wave (DDW) state [4], which was proposed 
Jit] as a candidate state responsible for the many anomalous 
properties of the pseudogap phase. Our results on diamag- 
netism will also be important in light of the recent experiments 
in Ref. [6] which point towards an alternative source different 
from vortices to explain the large diamagnetism observed in 
the pseudogap state of the cuprates. 

The ordered DDW state consists of counter propagating 
bond currents on the neighboring plaquettes of a 2D square 
lattice (Fig. 1) |4[, which can be taken as the Cu lattice of 
the high-T c cuprates [5]. The diamagnetism of this state has 
already been examined within a mean field description J3l 
in which the direction of the currents on the bonds remains 
frozen. The only source of diamagnetism in this descrip- 
tion are the nodal quasiparticles, whose contribution has been 
shown to be exceedingly small [7]. However, the mean field 
description does not include the direction fluctuations of the 
bond currents themselves. Because fluctuating bond currents 
respond much more strongly than quasiparticles to an applied 
magnetic field (see below), it is possible that these fluctuations 
gives rise to an enhanced diamagnetic response. A suitable 



way to include these direction fluctuations is to formulate the 
DDW state in terms of a vertex model, in which the directed 
arrows represent directed bond currents (Fig. 1). In this pa- 
per we use this vertex model description of the DDW state 
to show that the diamagnetism of the state significantly en- 
hances with increasing temperatures. Further, including also 
the magnitude fluctuations of the bond currents (not contained 
in the usual 6-vertex model), which are important at high tem- 
peratures, we show that the high-T, low-field, x can diverge 
as a power law of the applied field H. 

6-vertex model and F-model The classical vertex models 
were originally proposed to study anti-ferroelectric materials 
and associated phase transitions in electric fields UGH. One 
specific vertex model, called the 6-vertex model, is particu- 
larly interesting since it can be solved exactly by transfer ma- 
trices ifsL [Til fl2h . The 6-vertex model is defined by a set of 
vertices constructed out of directed arrow variables defined on 
the bonds of a square lattice. The arrows can represent any di- 
rected classical variable which serves as the building block of 
a thermodynamic statistical mechanical system. To describe 
an orbital current system the arrows are taken as directed bond 
currents. On the 2D x — y plane, each of the nearest-neighbor 
bonds in the 4 directions d = ±ax, ±ay from a vertex v is as- 
sociated with an orbital current Z„ ^ of magnitude Iq. The cur- 
rent Iv is positive for current flowing parallel to d and neg- 
ative for current flowing anti-parallel to d. In the steady state, 
there is no charge accumulation at each vertex and therefore 
the current is divergence-free (^2 d Iv = 0). This mandates 
that the total number of possible vertices on a square lattice is 
WjIj = 6 (Fig. 1). The Hamiltonian describing the orientation 
of the currents for a special case of the 6-vertex model called 
F model (0 U ])is given by 
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As is clear from Eq. Q] in the F-model the anti-ferroelectric 
(AF) vertices (defined by 1% = —Iy X and if = — I~ v ) are 
assigned negative energies —K and the rest of the vertices 
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FIG. 1. (a) The six possible current vertices in the 6- vertex model. 
The vertices (5) and (6) are the AF vertices energetically favored by 
the F model, (b) The AF aligned low-temperature ground state of the 
F model corresponds to the d— density wave (DDW) state proposed 
for the cuprate superconductors. 



have energy 0. Therefore, at low T, the ground state of the 
F-model is the ordered AF state, which is nothing but the or- 
dered DDW state when the arrows represent currents. The AF 
state survives thermal fluctuations up to a critical temperature 
T = T* . Above T* the current variables disorder into a criti- 
cal phase with power-law current correlations 10]. 

DDW state and its relevance to the high-T c cuprates. 
The singlet DDW state, described by an order parameter 

(fk+Q,a^,p) oc iW k 5 al 3, W k = -^(cosfcz - cosk y ), 

Ck,c k are fermion operators, k is a 2D momentum, Q = 
(tt, 7r), and a, f3 are spin indices, has been proposed as provid- 
ing a phenomenologically consistent explanation for the pseu- 
dogap phase of the underdoped cuprates [5]. The assumption 
of DDW order below optimal doping can lead to an explana- 
tion of numerous experiments including the abrupt suppres- 
sion of the superfluid density 11311 and Hall number [ 14] below 
optimal doping as well as the more recent quantum oscillation 
experiments [15] and Nernst effect [116]. Mathematically, any 
Hamiltonian that leads to d-wave superconductivity in the un- 
derdoped cuprates will almost certainly favor DDW order as 
well 1 4, 17], making their coexistence and competition in the 
phase diagram a plausible scenario. 

Connection of DDW state with F-model. In a mean field 
picture, the only way the DDW state can thermally disorder 
is via a collapse of the magnitude of the order parameter Wo 
(i.e. collapse of the magnitudes of the currents themselves) at 
a second order thermal phase transition. However, this mean 
field description does not take into account the possible direc- 
tion fluctuations of the bond currents. As is clear from Fig. 1, 
the ordered DDW state is nothing but the low-T AF state of 
the F model. In F model, with increase in T, the direction fluc- 
tuations of the currents eventually make the system pass into 




FIG. 2. (Color online) Two states, AF and the maximally current 
carrying state I of the 6-vertex model. The AF state is characterized 
by small clock-wise loops around plaquettes on one of the sublattices 
(marked by X). To construct the maximally current carrying state 
I, one starts with the AF state and reverses the counter-clockwise 
currents adjacent to all closed loops (red dashed curves). See text for 
details. 



a current disordered state above the temperature T*. Thus 
the DDW state can thermally disorder by bond current fluc- 
tuations above T* long before the order parameter magnitude 
Wo itself collapses at a mean-field temperature T m > T*. In 
this way the F-model and its quantum extension have recently 
been used HI 1911 to describe the thermal and quantum disor- 
dering transitions of the DDW state in the underdoped regime 
of the cuprates. 

Diamagnetic response of F-model. Let us first give an in- 
tuitive argument for the diamagnetism of the F-model. The 
interaction of the orbital currents with an external magnetic 
field can be described by a term H mag = — J drJ(r) ■ A(r) 
where J(r) is the orbital current density and A(r) is the vec- 
tor potential. The expression for H mag is derived by applying 
the minimal substitution p — > (p — qA/c) to the Schrodinger 
equation for the electrons. The divergenceless orbital cur- 
rent J(r) can be expressed in terms of a magnetization den- 
sity m(r), J(r) = V x m(r). Using the vector magnetic 
field B = V x A and integration by parts lead to the form 
H mag = - J drJ(r) ■ A(r) = J drm(r) ■ B. Thus, the 
bond currents can lower energy by aligning in circles perpen- 
dicular to B so that B ■ m(r) < (hence the response is 
diamagnetic). If such an arrangement of bond currents can 
be accessed by flipping the local vertices in closed loops (to 
maintain the charge conservation)(Fig. 2), the resultant dia- 
magnetic response can be large . Below we discuss this more 
quantitatively. 



Adding the magnetic interaction term H mag to the Hamil- 
tonian in Eq.[T]we get the total Hamiltonian as: 

H = £[-- (4 d) - It d) ) 2 alWd ■ A(v + |)]. (2) 

Here A(x,y) is the vector potential given by A(x,y) — 
^(—By, Bx). The second term is equivalent to H mag when 
the current densities are limited to the bonds and is equivalent 
to B ■ M where M is the total magnetic moment. To cal- 
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culate the magnetization density we write the partition func- 
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tion Z = £\ 

e~ E °^ kBT [l+Y^jjto e~ ( - E ^ E °^ kBT ], where Ej is the energy 
associated with the configuration j and j = is the minimum 
energy configuration. Using Eq. the energy Ej is given by 
_Ej = EJ + B Mj where Ej is the energy of the first term of 
H and Mj is the total magnetic moment of configuration j. 
Calling the configuration with the maximum diamagnetic mo- 
ment I (Fig. [2]right panel), the free-energy F — —ksTlogZ 
can be written as 



F = B Mi + E Q — k B T\og 



(3) 



where ^ = E? + B ■ [Mj - M{\ and Eq = mhijEj. Since 
E < Ei = E" and B ■ M > B ■ Mi, it follows that Eft < 
Eq < Ef . Additionally, the magnitude of the energy Ej of 
any state j from Eq.Q]must be less than KI 2 R 2 where R is 
the radius of the system containing ~ R 2 vertices. Thus the 
second term Eq is bounded by \Eq\ < KIqR 2 . The third term 
in Eq. [3]which is logarithmic also scales as R 2 since each of 
the terms under the summation over j is less than unity and 

r-,2 

there are at most 6 such terms corresponding to the current 
configurations on R 2 vertices. Combining these results, we 
find that 



F = BMi + 0(R 2 ) 



(4) 



where 0(R 2 ) represents corrections of order R 2 (which can 
be neglected as \M{\ ~ IoR 3 /aas we show below). 

The state / with maximum diamagnetic moment is under- 
stood by starting with the AF state as follows (Fig. [2]): Imag- 
ine large closed loops (red dashed curves) passing through the 
dual lattice points marked by the crosses. The currents on 
the bonds touching these loops but on the two opposite sides 
flow in opposite directions. For example, on the left panel of 
Fig. 12 the bond currents right outside the loops are clockwise 
and right inside are counter-clockwise. Reversing the counter- 
clockwise currents touching all such closed loops leads to the 
state I shown on the right panel. The total magnetization M is 
the product of the current and the total area enclosed by all the 
clockwise loops. Since most of the clockwise loops (~ R/a 
in number) enclose an area of order R 2 , the total magnetic 
moment is Mi oc —zI^R 3 /a which is the desired result. Us- 
ing this equation for Mi and neglecting terms of 0(R 2 ) for 
large R we get from Eqn.H] F(B) = —BIqR 3 /a. The mag- 
netization density is calculated as m(B) = -j^Qpg- This gives 
m(B) ~ —zI R/a, which is divergent in the thermodynamic 
limit (R — » oo) for any non- vanishing B. 

So far we have ignored the magnetic field generated by the 
induced currents themselves. Such a field results in a magne- 
tostatic current-current interaction. The current-current inter- 
action can be accounted for by using simple magnetostatics, 



H = B + Air\m(B) 



(5) 



where B is the magnetic field and H is the magnetic induction 
which can be taken as the externally applied magnetic field. 
Here m{B) is the magnetization density which is opposite in 
direction to B (diamagnetic) in sign and increases in magni- 
tude from \m(B)\ = at B = to \m(B)\ ~ IqR/cl for any 
magnetic field B > 0(1/ R). To estimate the solution B of 
Eq.|5]let us define the function f(B) = B + ATt\m(B)\ - H 
such that Eq. [5] is written as f(B) — 0. Since H < I^R/a 
(R — > oo in the thermodynamic limit), it follows that f(B) > 
OforS > 0(1/ R). On the other hand / (B = 0) = -H < 0. 
Therefore f(B), being an increasing function of B, has a 
unique root satisfying the constraint B < 0(1/ R). Since 
in the thermodynamic limit R — > oo, it follows that B is com- 
pletely expelled from the system and it behaves like a perfect 
diamagnet similar to a type I superconductor. 

Diamagnetism in the AF phase. Despite the above analy- 
sis, the low-T (T <C T*) AF phase being gapped is not ex- 
pected to have a large diamagnetic response. The response 
of the AF phase to a magnetic field should be dominated by 
flips of small current loops. The combination of these ele- 
mentary current loop flips can generate a flip of a large loop 
of length L which has a magnetic moment \ M\ ~ IqL 2 . From 
Eq.Q]such loop flips cost energy KL (flipping each AF vertex 
costs energy K and L such vertices need to be flipped for a 
loop of length L). However, the applied field lowers the en- 
ergy of such a current loop by —B\M\ = —BI$L 2 . Thus 
the energy cost of a flipped loop, V(L) = KL — BIqL 2 , is 
positive for small L, has a positive peak at L = K/2BIq 
and becomes negative for large L. Thus only loops that 
form out of a thermal fluctuation with an energy larger than 
maxLV(L) — K 2 /2BIq can cross the threshold value of 
L = K/2BIq to become a large loop. The fraction of such 
high energy loops is determined by the Boltzmann factor as 

e ~K 2 /2k B TBI ^ Thus the A p state is stab j e for l ow .rp and B 

(i.e e~ K l 2k BTBi ^ y) This conclusion is also consistent 
with the numerical monte-carlo simulations [23]. It also fol- 
lows that as T increases (T is a substantial fraction of T*), the 
diamagnetic response of the AF state should increase signif- 
icantly. However, the time-scale for development of diamag- 
netic response is expected to become longer as temperature 
becomes smaller than T* leading to possible hysteretic be- 
havior of the magnetization as a function of applied magnetic 
field. 

Diamagnetism in the critical phase. For T > T* , the AF 
order is completely destroyed and the system develops crit- 
ical current fluctuations whose correlation is scale-invariant. 
From the argument in the previous paragraph it follows that 
this phase is strongly diamagnetic. The scale invariance of 
the fluctuations in the critical phase allows one to describe the 
critical phase by a continuum theory such as a height model 
I 2l[|22ll . In the height model the magnetization density m(r) 
is mapped to the vertical displacement h(r) of a 2D surface 
such that m(r) = Ioh(r)z and J(r) = V x m(r). 

The current-current correlations in the high T critical phase 
of the F model are obtained from the Gaussian theory of height 
fluctuations 12 ll 12211 described by the coarse-grained contin- 



4 



uum Hamiltonian 

H= [ d 2 rK\\7h\ 2 — B ■ M u 



(6) 



where Mtot as before is the total magnetization of the model. 

From the argument in the last section, it is clear that for 
T > T* , the F model responds to a magnetic field by gen- 
erating large current loops. This results in the formation of 
patches of circulating currents. If the sizes of these patches 
are macroscopically large, this would lead to perfect diamag- 
netism. However, so far we have neglected the magnitude 
fluctuations of the bond currents Iq, which should be taken 
into account at high T. Such magnitude fluctuations of Iq can 
occur from spontaneous thermal fluctuations of the DDW gap 
magnitude Wo and variations of the local density of quasipar- 
ticles. Here by quasiparticles we mean the quasiparticles in 
the DDW state that carry charge ll24ll . In the presence of such 
magnitude fluctuations the bond-current Iq is allowed to vary 
spatially as Io(r,t) so that the current density is now given 
by J(r, t) = Io(r, t) V X (h(r)z). Therefore it is no longer 
a divergence-free quantity. The magnitude fluctuations of the 
currents will lead to a cut-off length-scale Rq > R P for the 
patch sizes R p . The introduction of this cut-off length-scale 
in the original F model directly leads to a power-law depen- 
dence of m on B as we show at the end of this section. Below 
we first speculate on a possible mechanism for the emergence 
of this cut-off scale. 

Let us consider a patch with a circularly symmetric cur- 
rent density J(r,t) = I (r,t)(Vh(r) xz) = l o (r,t)h'(r)0 
corresponding to a circularly symmetric height profile h{r). 
Apart from a thermal fluctuation component, £(r, t), the mag- 
nitude of the local current Io(r, t) depends on the density of 
quasiparticles n(r,t). Writing Io(r,t) = Io(S(n(r,t)) + 
C(r, t)), the time-varying current density J(r, t) is now given 
by J{r,t) = l o (r,t)h'(r)0 = I h'(r)[S(n(r,t) + ■ 
V({r,t)]0 = J{r)[S{n(r,t))+0-V({r,t)}0, where h'(r) = 
dh J^ , C( r iO i s tr, e noise term that accounts for the sponta- 
neous thermal fluctuations of the current density and is the 
tangential direction around a loop on the circular patch. Since 
charge density is conserved we have the continuity equation, 
dtn(r, t) + V • J(n(r, t)) — D\7 2 n, where D is the diffusion 
constant of the quasiparticles. For a circular profile of a single 
patch of circulating currents, this equation can be written as 

d t n + J(r)[S' (n)0 ■ Vn + ■ V((r,t)} = DV 2 n. (7) 

The second term in the above equation, which is referred to 
as the bond-current term drives current only along the tangen- 
tial direction to the loop at a given radius and hence is one- 
dimensional in character. The quasiparticle diffusion term 
on the right hand side of the above equation is two dimen- 
sional in space and in general will have a finite radial com- 
ponent. However, long wave-length fluctuations in the bond 
current £(r, t) (the correlation function of £(r, t) is taken as 
(C(r, t)C(r', t')) ~ S{r - r')5(t - t% with a length-scale A, 
create long wave-length variations in the quasiparticle density 



which are tangential in direction (as shown by the third term 
in Eq. |7). The evolution of the quasiparticle density follow- 
ing such a fluctuation, which is described by Eq. [7j is domi- 
nated by the tangential bond-current term (i.e. second term), 
which scales as A -1 , and the contribution of the radial dif- 
fusion term (which scales as A~ 2 ) is subdominant for large 
A. Therefore, the long length-scale behavior of Eq. [7] can be 
understood in terms of approximately decoupled ID conser- 
vation law equations for each loop at radii \r\ = r. As shown 
m Refs. E E 12811 . for mean-field bond-current densities 
S(n) which vanish for very small or large quasiparticle den- 
sities, such ID equations are unstable to the proliferation of 
long-lived shock solutions and the current density Iq is ex- 
pected to drop to zero at some radius Rq which can be de- 
termined by solving the relation J(r — Rq) cx Rq". In 
the present case, how such a relation arises can be seen in 
the following way. From Eq. [71 we see that the bond-current 
term at a length-scale Rq scales as J(r = Rq)Rq 1 , while 



the diffusion term scales as DR n 



Thus the bond-current 



term, which drives the instability, dominates at a radius Rq 
when J(r = Rq)R q ~ 1 ~ DRq 2 . This leads to the constraint 
J(r = Rq) cx Rq° with a = 1. In reality this simple argu- 



ment ignores the Rq -dependence of n s 125L 12711 . which leads 
to a value of a slightly greater than 1 . 

Let us now return to the F model and show that the intro- 
duction of a cut-off scale leads to a power-law dependence of 
to on B. We consider the free-energy of a patch of finite ra- 
dius Rq. Using Eq.|6]the free-energy density of such a patch 
is given by 

" Ro 2rdr , ~ 



f 



■{Kb! 2 -IqBK). 



(8) 



Minimizing / with respect to h using the usual method of 
variations (and without any approximations) leads to h'(r) — 
-BIqt/AK. Using the relation I h'(Ro) = Rq", we get 
R cx Substituting this in Eq. [8] leads to / cx 

B 2 Rl - B 2a ^ 1+a \ The diamagnetic susceptibility x = 
cx _B~ 2 /( 1+a ) diverges as an inverse power-law. 

Conclusion. We analyze the diamagnetic response of the 
6-vertex model, which is used to model the DDW phase pro- 
posed for the pseudo-gap phase of the high-T c cuprates. We 
find that deep in the low-T AF phase the diamagnetic response 
is weak. With increasing T, especially for T <T* , the dia- 
magnetism is significantly enhanced. The disordered critical 
phase for T > T*, is perfectly diamagnetic within the strict 
6-vertex model. With the magnitude fluctuations of the cur- 
rent (magnitude fluctuations of the DDW order parameter Wq) 
taken into account, the low-field diamagnetic susceptibility x 
in this phase diverges as a power-law of the field. 
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